Quantum Monte Carlo and exact diagonalization study of a dynamic Hubbard model 
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A one-dimensional model of electrons locally coupled to spin- 1/2 degrees of freedom is studied by 
numerical techniques. The model is one in the class of dynamic Hubbard models that describe the 
relaxation of an atomic orbital upon double electron occupancy due to electron-electron interactions. 
We study the parameter regime where pairing occurs in this model by exact diagonalization of small 
clusters. World line quantum Monte Carlo simulations support the results of exact diagonalization 
for larger systems and show that kinetic energy is lowered when pairing occurs. The qualitative 
physics of this model and others in its class, obtained through approximate analytic calculations, is 
that superconductivity occurs through hole undressing even in parameter regimes where the effective 
on-site interaction is strongly repulsive. Our numerical results confirm the expected qualitative 
behavior, and show that pairing will occur in a substantially larger parameter regime than predicted 
by the approximate low energy effective Hamiltonian. 



I. INTRODUCTION 

Dynamic Hubbard models have been recently intro- 
duced as a new class of model Hamiltonians to de- 
scribe the relaxation of atomic orbitals when electrons 
are added to orbitals already occupied by other electrons 
. This process, originating in the strong on-site re- 
pulsion between electrons in the same atomic orbital, is 
not described by the conventional Hubbard model |Q. 
In dynamic Hubbard models this physics is represented 
either by introducing auxiliary spin || or oscillator ||] 
degrees of freedom, or by adding a second electronic or- 
bital to the site Hilbert space 0, with suitable interac- 
tion parameters. As a consequence, the on-site Hubbard 
repulsion becomes a dynamical variable and can take a 
range of values rather than a single fixed value as in 
the static (conventional) Hubbard model. It has been 
proposed that this physics is ubiquitous to electrons in 
atoms, molecules and solids and that it is rele- 

vant to the understanding of superconductivity in nature 

While a vast amount of work has been performed over 
the years on the conventional Hubbard model fnio , 
very little work has been done so far on dynamic Hubbard 
models. It is known 12 (^JtJ that in the strong coupling 
anti-adiabatic limit these models map onto the Hubbard 
model with correlated hopping, i.e. a Hubbard model 
where the electronic hopping amplitude depends on the 
occupation of the two sites involved in the hopping pro- 
cess. This model is known to exhibit superconductivity 
when the Fermi level is close to the top of the band, both 
from mean field calculations [fU3|-[l5|] , exact diagonaliza- 
tion , and other exact techniques [|l9|j2C| 1 . Further- 
more, a variety of observable properties have been calcu- 
lated in this limit such as thermodynamics |fl3| , ^l|| , tun- 
neling p2[, optical properties pa], pressure dependence 



|2l| , etc. Because superconductivity occurs in the dilute 
carrier concentration regime it is believed that these BCS 
mean field calculations are reliable pi| . 

The antiadiabatic limit of these models occurs when 
the frequency of the associated boson degree of freedom, 
ljq , is much larger than the effective hopping amplitude 
for the electrons (small polaron regime (25)), where the 
boson follows the electronic motion. In that limit the 
parameter regime where pairing occurs can be calcu- 
lated exactly for a dilute concentration of hole carriers 
p4j . Furthermore, numerical calculations on finite clus- 
ters show that the doping regime where pairing occurs 
is accurately estimated by BCS theory |l8|. For finite 
frequency loq, some numerical results have been reported 
P,^2|. However it is generally not known whether finite 
ujq enhances or reduces the tendency to pairing. 

Furthermore, in the antiadiabatic limit the single car- 
riers have large effective mass, and the effective mass is 
lowered when carriers pair ^g,^7|. The resulting gain 
in kinetic energy drives superconductivity ]28[ ]. It is not 
known whether this physics exists beyond the antiadia- 
batic limit. 

In this paper we study a particular realization of a 
dynamic Hubbard model, with an auxiliary spin degree 
of freedom, by exact diagonalization of small clusters 
and a quantum Monte Carlo method, to shed light onto 
the properties of the model away from the antiadiabatic 
limit. We believe that similar qualitative behavior may 
be found in the entire class of dynamic Hubbard mod- 
els. Briefly, our results show that the qualitative physics 
of the antiadiabatic limit persists for finite loq, and that 
the parameter regime where pairing occurs can be sub- 
stantially larger. Even though our results are for a one- 
dimensional system, we believe it is likely that the same 
occurs in higher dimension. 

The model studied here bears some superficial resem- 
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blance to electron-boson models that have been exten- 
sively studied in the past such as the Holstein model 
p5| . However its physics is qualitatively different. To 
illustrate this point we present some numerical results 
for an electron-hole symmetric model with an auxiliary 
spin degree of freedom coupled to the electronic site den- 
sity. This model is expected to be similar to the Holstein 
model, and exhibits qualitatively different physics to the 
dynamic Hubbard model. 

The paper is organized as follows. Section II defines 
the model and discusses its properties in the antiadia- 
batic limit. In Sect. Ill we present results for the effec- 
tive interaction and kinetic energy from diagonalization 
of small clusters, and Sect. IV discusses results of world 
line quantum Monte Carlo simulations. In Sect. V we 
present and discuss results for the electron-hole symmet- 
ric Holstein-like model. Sect. VI discusses the relation 
between the dynamic Hubbard model studied here for a 
site and a real atom. We conclude in Sect. VII with 
a summary of our results and a discussion of the many 
open questions in this area. 



A. Site Hamiltonian 

We will use the Hamiltonian in the hole representation, 
Eq. (2). The site eigenstates when there are n holes at 
the site are, in terms of the spin-1/2 a z eigenstates |+ >, 
|-> 



\n >= u(n)\+ > +v(n)\- > 



\n >= v(n)\+ > -u(n)\— > 



(4a) 



(4b) 



with eigenvalues (excluding the er-independent term in 
Hi) 



and 



5 (n) = -e(n) = ~uj y/l+g 2 



« 2 (0) = tt(1 + 



2 v/T+? 



(5) 



(6a) 



II. DYNAMIC HUBBARD MODEL WITH 
SPIN-1/2 DEGREE OF FREEDOM 



^ 2 (0) = -(i 



(6b) 



The essence of dynamic Hubbard models is electron- 
hole symmetry breaking at the local (1-site) level, so that 
the dressing of a hole is larger than the dressing of an 
electron || . This physics originates in the dynamic low- 
ering of the on-site repulsion U when a second electron is 
added to an atomic orbital, due to rearrangement of the 
first electron, and is a ubiquitous phenomenon in atoms 
|^| . There are a variety of dynamic Hubbard models that 
can be constructed with an auxiliary spin-1/2 degree of 
freedom . Here we consider the site Hamiltonian for 
electrons 

Hi = lj <t x + guj Q a l z + [U - 2guj al}n^n ll (la) 

Hence, for zero and one electrons (one and two holes) at 
the site, the site Hamiltonian is 



u(0>(0) 



Hi(rii < 2) = LO Q al + guj a l z 



(lb) 



and for two electrons (zero holes) at the site it is (spin 
part only) 

Hi{rii = 2) = uj a x - gu a a\ (lc) 

Eq. (la) can be written in hole representation as 

Hi = u) a x + goj [2(n^ + mi) - 

+[U - 2guoo- i z ]n n n il (2) 

(omitting a chemical potential term), and the lattice 
Hamiltonian is 



H = J2H i -tJ2[clc i+1 , l7 +h.c.} 



(3) 



in either electron or hole representation. The electron- 



and 



2^/^TF 



t (0) 



= 9- V^+g 2 



v(0) 

u(l) = tt(2) = u(0) 
= v (2) = u(0). 



(6c) 



(6d) 



(7a) 



(7b) 



hole transformation is 



Hence the ground state energy is independent of the elec- 
tronic site occupation in this model. The site cigcnfunc- 
tions depend on g but not on u>q, and are the same for 
site occupation n — 1 and n = 2, and different for n — 0. 
For large g the ground state wavefunctions are almost 
eigenstates of a z , with a z ~ — 1 for one-hole and two- 
hole occupation, very different from the one for zero hole 
occupation for which a z ~ +1; while for small g the 
ground state wavefunction is almost an eigenstate of a x 
{o~x ~ — 1) and similar for the different hole occupations. 
The site eigenvalues depend on both g and c^o- 

The on-site repulsion between two holes (or two elec- 
trons) at the same site depends on the state of the spin 
degree of freedom and can range between U + 2gujQ and 
U — 2gu>Q. The effective on-site repulsion however, since 
the ground state energy Eq. (5) is independent of oc- 
cupation, is simply U . Note that our notation here is 
different from that of Ref . j2| , where U denoted the bare 
on-site repulsion; here, the bare on-site repulsion, which 
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is the on-site repulsion if the background degree of free- 
dom is not allowed to relax upon double occupancy, is 



U ba 



u - 



2g 2 u 



(8) 



or Ubare ~ U + 2gujQ for large g. Finally, the overlap 
matrix elements between the ground state wavefunctions 
for the various hole occupations are 



< 0|1 >= 2u(0)u(0) 



1 



< 112 >= 1 



(9a) 



(9b) 



B. Effective low energy Hamiltonian 

The effective hopping amplitude for a hole between 
neighboring sites if the spin degree of freedom makes a 
ground-state to ground-state (diagonal) transition is 



to 



< o|i > \ 2 t = SH 



t 



i 



(10a) 



if there arc no other holes in the two sites involved in the 
hopping process. Instead, if there are either one or two 
other holes of opposite spins the hopping amplitudes are 



h = | < 0|1 >< 1|2 > | = St 



t = I < 112 > rt = t 



(10b) 



(10c) 



respectively. The low energy effective Hamiltonian for 
holes in the small polaron regime is then 

2,0" i 

(11a) 

t% = t[S 2 + S(l - S)(n^ a + nj ,- a ) + (1 - Sfn^nj. _ c 

(lib) 

In the regime of low hole concentration the hopping pro- 
cesses where more than two holes are in the sites involved 
can be neglected, and the effective Hamiltonian is the 
Hubbard model with correlated hopping |l3| 

H e ff = - + ^■ t ( n t,-<y + ™;+i,-<x)K c L c i+i<x + h.c.) 

+U^2n^mi (12a) 



At = h -t 2 =tS{l- S) 



(12b) 



The binding energy of the polaron is obtained from the 
difference of e(0) and the expectation value of Hi if the 
spin doesn't adjust to the presence of a carrier, and yields 



2uJ g 2 



,9 Z 



(13) 



The criterion for small polaron formation is that the en- 
ergy of the polaron is smaller than that of a carrier that 
moves without changing the spin background: 



z(t - t 2 ) < e p 

with z the number of nearest neighbors to a site (z 
in one dimension). Eq. (14) yields 



,9 Z 



(14) 
= 2 

(15) 



as the condition for polaron formation. For the hopping 
of a single polaron, the antiadiabatic limit is valid if the 
polaron hopping amplitude t 2 is smaller than the spacing 
between site energy levels, hence 



t< Wo (l+5 2 ) 3/2 



(16) 



which is always satisfied if the condition Eq. (15) is sat- 
isfied. 

The condition Eq. (15) indicates that the small po- 
laron regime will occur when either w is large or the 
coupling g is large. However, the condition Eq. (15) is 
not sufficient for the effective Hamiltonian Eq. (11) or 
(12) to be accurate in the presence of more than one car- 
rier. Virtual transitions of a hole to a nearest neighbor 
site occupied by another hole yield a contribution to the 
effective interaction between holes with an amplitude of 
the form 



< 110 > 



1 



5(0) + e(2) - 2e(l) + U 1 + g 2 U + 2lo ^/TTg~ 2 



(17) 



from 'vertical' transitions, which can be much larger than 
the second order contribution from the effective Hamil- 
tonian Eq. (11) that describes only diagonal transitions: 



< 110 > 



1 



e(0) + e(2) - 2c(l) + U (1 + g 2 )U 



(18) 



These contributions from vertical transitions can be ne- 
glected if Eq. (17) is smaller than Eq. (18), which yields 
the condition 



w 



y/1 + g 2 > U 



(.9 2 - 1) 



or, for large g 



w > 



Ug 



(19a) 



(19b) 



Only when both conditions Eq. (15) and Eq. (19) hold 
can the effective low energy Hamiltonian Eq. (11) or 
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(12) be expected to be accurate. In particular, for large 
g and small ojq Eq. (15) may hold and Eq. (19) may not. 
In that case, the effective Hamiltonian Eq. (11) can be 
expected to underestimate the tendency to pairing due 
to its complete neglect of the site excited states. The 
antiadiabatic limit where the effective Hamiltonian Eq. 
(11) is valid hence occurs for ujq — > oo for fixed g but not 
for g — > oo for fixed u>q. As g increases it is seen from 
the condition Eq. (19) that the antiadiabatic limit will 
be attained for larger lvq. 



C. Pairing condition and effective mass in 
antiadiabatic limit 

The condition on the parameters of the Hamiltonian 
Eq. (12) to yield pairing of two holes in a full band is (in 
one and two dimensions) p3] 



At 



> \1 + 



U 



1 



(20) 



with L\ = 2zt 2 the single carrier renormalized band- 
width. This is also the condition for superconductivity 
within BCS theory in the dilute limit in any dimension 
[n~3|. Using Eqs. (12b) and (10a) it translates into 



U 



< 



9 



(21) 



which shows that for g > 1 pairing will occur even if the 
on-site repulsion is larger than the effective bandwidth. 
Eq. (21) can also be written as 



U g 2 
D ~ 1 + 



(22) 



with D = 2zt the unrcnormalized bandwidth. 

The polaron hopping amplitude increases as the hole 
filling of the band increases, according to 



t(nh) = t 2 + ntAt 



(23) 



with rih the average number of holes per site (0 < < 
2) , and correspondingly the bandwidth increases 



D{n h ) = D h {l + nA 

t2 



(24) 



from D(n h = 0) = D h to D(n h = 2) = D. The polaron 
effective mass correspondingly decreases as the number 
of holes increases 



m*(n h ) = 



2t(n h )a 2 



(25) 



with a the lattice spacing. 

When two holes bind in a pair, the pair hopping ampli- 
tude t p in the dilute hole concentration regime is found 



to be always larger than 1/2 the single particle hopping 
amplitude 



t P > t 2 /2 



(26) 



i.e. the pair effective mass is smaller than the sum of the 
effective masses of its constituents [^7| . This is opposite 
to what happens in other models such as the attractive 
Hubbard model. Expressions for the pair mobility t p are 
given in refs [ p6| , p7| . The pair mobility is defined in terms 
of the energy dispersion relation for a pair of center of 
mass momentum q 



E(q) = E + t p q 2 



(27) 



and can be obtained by calculating the London penetra- 
tion depth in the dilute limit. The kinetic energy per two 
holes in the dilute limit when there is no pairing is 

< T s >= -4i 2 (28a) 

and when there is pairing the kinetic energy per pair is 



<T P >-- 



for a one-dimensional chain. 



-8t 



(28b) 



III. EXACT DIAGONALIZATION RESULTS 

The Hamiltonian of interest has 8 states per site, so 
that clusters of up to 8 sites could be studied with cur- 
rent computer capabilities. In this initial study we re- 
strict ourselves to 2 and 4 sites only. The results are 
qualitatively similar and we expect similar qualitative 
results for larger clusters, although quantititative differ- 
ences may be expected for weak coupling. We compute 
the effective interaction for two holes in a cluster from 
the usual formula 



U eff =2E (l)-E (p)-E (2) 



(29) 



with Eo(rih) the ground state energy for rih holes; U e ff < 
signals a tendency to pairing and superconductivity. 
For the N — 2 cluster the effective interaction in the 
antiadiabatic limit is 



U, 



eff 



ti = t 2 + At = tS 

and the condition for pairing (U e ff < 0) is 

U g 2 
— < — — 

2t ~ 1 + q 2 



(30a) 



(30b) 



(31) 



so that in the antiadiabatic limit pairing cannot occur 
for U > 2t for any value of the coupling parameter g. 
Throughout this and the following section we will use 
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units so that t = 1. For the N = 4 system, and in fact for 
any N > 4 the condition for pairing in the antiadiabatic 
limit is 



U g 2 
At~ 1 + 



(32) 



so that pairing will not occur for U > 4 for any g in the 
limit luq — > oo. 



A. Results for effective interaction 

Figure 1 shows the effective interaction for the N = 4 
cluster as function of coupling constant g, for various val- 
ues of the on-site repulsion U and two values of the fre- 
quency ujq, together with the results in the antiadiabatic 
limit. Note that for small (Figure la) the effective 
interaction is substantially more attractive than in the 
antiadiabatic limit. As ojo increases (Fig. lb) the results 
approach those of the antiadiabatic limit, as expected. 
The behavior of U e f f versus g is non-monotonic particu- 
larly for small values of U. 

In Figure 2 we show the dependence of the effective 
interaction on u>q for the N = 4 cluster for fixed U = 4 
and various values of g (a), and for fixed g for various 
values of U (b). The limiting values for ojq — ► oo are also 
shown (dashed lines). For [7 = 4 there is no pairing in 
the antiadiabatic limit, while Fig. 2a shows that for finite 
frequency pairing will occur for g > 2. Similarly, Fig. 2b 
shows that for g = 3 (corresponding to an effective mass 
enhancement m* /m = 1 + g 2 — 10) pairing will occur 
up to at least U = 6 at finite ujq, while in the antiadi- 
abatic limit U = 3.6 is the maximum on-site repulsion 
that allows pairing for g — 3 according to Eq. (32). Note 
that for larger g the antiadiabatic limit is approached for 
larger ujq, in accordance with the discussion following Eq. 
(19). 

Even a cluster as small as TV = 2 shows behavior rep- 
resentative of larger clusters and of (we believe) the ther- 
modynamic limit. The reader can easily verify that the 
effective interaction for the N = 2 cluster obtain by ex- 
act diagonalization closely resembles the behavior of the 
4-site cluster shown in Figs. 1 and 2. In the antiadi- 
abatic limit the effective interaction as function of g is 
monotononically decreasing with g if the condition for 
pairing Eq. (31) or (32) is not satisfied, while if it is 
satisified it has a (negative) minimum for a finite g that 
decreases as U decreases below the limits given by Eqs. 
(31) and (32). 

For an infinite chain, the pair binding e nerg y can be 
calculated exactly in the antiadiabatic limit [gTJ . The ap- 
pendix of Ref. 27 gives an analytic expression for the pair 
binding energy eb in one dimension. The quantity U e f t 
defined by Eq. (29) calculated here should go to — e& as 
the cluster size increases. Figure 1 also shows results for 
—eb (dotted lines), which go to zero when the parameters 
satisfy the equality in condition Eq. (32). The difference 



between the dotted and dashed lines gives the magnitude 
of finite size effects for the N = 4 cluster. It can be seen 
that the qualitative behavior of — ej, for the infinite chain 
and U e f f for the 4-site chain is the same. The effect of fi- 
nite size is to give a somewhat larger attraction, however 
the condition for pair formation {Ueff < 0) is the same 
for the N=4 cluster and the infinite chain in the antiadi- 
abatic limit (the dashed and dotted lines in Figure 1 go 
to zero at the same value of g). 

In Figure 3 we show the phase diagrams for the N = 2 
and N = 4 clusters indicating the region where pairing 
will occur for some finite frequency in this model. The 
full lines show the results in the antiadiabatic limit, Eqs. 
(31) and (32). It can be seen that the region of parameter 
space where pairing occurs is substantially enlarged for fi- 
nite frequency. Figure 4 shows the optimal frequency for 
pairing at the phase boundary for pairing , for the N = 2 
and N = 4 clusters. For large g, ujq increases slowly with 
g and is between 1 and 2 (in units of t). For decreasing g, 
ujq goes through a minimum and then diverges, in accor- 
dance with the fact that the phase boundary lines in Fig. 
3 merge with the ones in the antiadiabatic limit as g — > 0. 
It should also be noted that for points away from the 
phase boundary the optimal frequency that gives maxi- 
mum attraction can be considerable smaller than those 
shown in Fig. 4 (see e.g. Fig. 2a for U = 4, g = 4 or 
Fig. 2b for g — 3, U — 2 where the optimal frequency is 
cjo-0.5). 



B. Results for kinetic energy 

The condensation energy in this model is known to be 
provided by lowering of kinetic energy in the antiadia- 
batic limit. Exact expressions for the pair kinetic energy 
for two bound holes in a one-dimensional chain are given 
in Ref. ^7|. Figure 5a shows exact results for the pair 
kinetic energy versus coupling constant g for the effective 
Hamiltonian Eq. (12) for various values of U. Pairing oc- 
curs for couplings obeying the condition Eq. (32), which 
for U — 0.8, 2, and 3.2 corresponds to g = 0.5, 1 and 2 
respectively. For g larger than those values the kinetic 
energy is given by the dashed line, lower than the full line 
which would be the kinetic energy in the absence of pair- 
ing. Note that even though the kinetic energy of a pair is 
lower than that of the unbound holes, it still decreases 
in magnitude as the couping g increases. Instead, the 
kinetic energy lowering, i.e. the difference between the 
kinetic energy of the pair and of the unbound holes, is 
non-monotonic, peaking at an intermediate g, similarly 
to the pair binding energy given by —U e ff. Figure 5b 
shows the kinetic energy lowering per pair 



AT 



--<T p > 



<T S > 



(33) 



and U e ff for the infinite chain, which is the negative of 
the pair binding energy eb calculated in Ref. |27|. It can 
be seen that the two quantities follow similar behavior 
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with coupling. In fact, their ratio is essentially constant 
as function of g for large g, as shown in Fig. 5c. As 
the pair binding decreases, either because U increases or 
g decreases, both U e ff and the kinetic energy lowering 
go to zero. However U e ff approaches zero quadratically 
|24f while the kinetic energy lowering approaches zero 
linearly, hence their ratio diverges as the pair binding 
energy goes to zero. 

Note that the kinetic energy lowering upon pairing is 
always lower than — e;,. This indicates that the potential 
energy change is positive, that is, there is a potential 
energy cost upon pairing, given by 

AU pot = (< T s > - < T p >) - e b (34) 

and the pair binding energy is smaller than would be ex- 
pected from the magnitude of kinetic energy lowering. 
The potential energy cost arises from the increased ef- 
fect of the on-site repulsion between members of a pair 
since the pair wavefunction has higher probability for site 
double occupancy. 

In the infinite chain there is a sharp phase transition 
between the state where the pair is bound and where it 
is unbound, indicated by the points in Figure 5a where 
the dashed lines join the full line. In the finite chain of 
course there is no sharp transition but rather a smooth 
crossover. Figure 6a shows results for the kinetic energy 
for a pair of holes in the 4-site chain in the antiadiabatic 
limit compared to the results for the infinite chain. As the 
coupling constant increases the 4-site results cross over 
from the kinetic energy of unbound holes to the kinetic 
energy of the paired holes. When g goes to zero the 
kinetic energy of two holes in the 4-site chain is slightly 
higher than the one for two unbound holes because of 
the effect of the on-site repulsion U ; this is of course 
a finite size effect, and for larger clusters and a fixed 
number of holes it will become negligible in the regime 
where the holes are not bound. Figure 6b compares the 
kinetic energy lowering and the pair binding energy for 
the 4-site chain in the antiadiabatic limit and the infinite 
chain for one case; it can be seen that both quantities 
follow similar behavior, and both are larger in magnitude 
than for the infinite chain. Note that in the finite chain 
kinetic energy lowering goes to zero for a smaller g than 
where U e tf goes to zero; this would also occur in the 
infinite chain for finite hole density. We conclude from 
these results that the 4-site cluster is appropriate to learn 
about the qualitative behavior of the kinetic energy for 
the infinite chain just as well as it is for the pair binding 
energy. 

Hence we can now learn about the effect of finite fre- 
quency on kinetic energy lowering by studying the 4-site 
chain. Figure 7a shows results for a finite small frequency, 
luo = 0.5, compared to the antiadiabatic limit luq = oo. 
Just as for the effective interaction(Figure la), the ki- 
netic energy lowering can be substantially larger for fi- 
nite frequency than for infinite frequency, and the largest 
kinetic energy lowering occurs for larger g for small fre- 



quency. Similarly Figure 7b shows kinetic energy lower- 
ing as function of frequency for fixed g. Similarly as the 
corresponding results for effective interaction Figure 2b, 
the kinetic energy lowering is largest in magnitude at a 
fairly low frequency and is considerably larger than in the 
antiadiabatic limit. The ratio of kinetic energy lowering 
to effective interaction for finite frequencies behaves sim- 
ilarly as in the antiadiabatic limit; this is shown for one 
case in Figure 7c. 

In addition to kinetic energy lowering it is of interest 
to consider the effect of finite loq on kinetic energy itself. 
Figure 8 shows results for single hole and pair kinetic en- 
ergy for finite ujq as function of the on-site repulsion U, 
compared with the limiting case ujq — oo. The difference 
between the dashed and full lines is the kinetic energy 
lowering. It can be seen that relatively speaking the ki- 
netic energy lowering is largest in the antiadiabatic limit, 
even though it is larger in magnitude for finite ujq. 

In summary, we have seen that the effect of finite fre- 
quency is to enhance the pair binding energy and the 
kinetic energy lowering found in the antiadiabatic limit. 
From these results we conclude that the pair condensa- 
tion energy in this dynamic Hubbard model also origi- 
nates in kinetic energy lowering, i.e. 'undressing'. As 
implied by the conductivity sum rule |^,^0|, lowering 
of kinetic energy should be accompanied by transfer of 
optical spectral weight from high to low frequencies, as 
well as by transfer of spectral weight in the single particle 
spectral function from high to low frequencies (3l]]. The 
one and two-particle spectral functions for the model Eq. 
(3) will be discussed in a separate paper. 

IV. QUANTUM MONTE CARLO SIMULATIONS 

With quantum Monte Carlo (qmc) methods one can 
study much larger systems than with exact diagonaliza- 
tion. We use the basis of er z eigenstates for the spin 
degrees of freedom, so that at every time slice i there are 
classical spins (Jj{i) at every lattice site j. The partition 
function is 

Z = Tre-? H = Tr]]e- ATH = 

i=l 

Tr ti E e- A ™^^ (35) 

i=l a j (i)=+/-l 

with At — (3/L, L the number of time slices. There are 
two basic approaches to quantum Monte Carlo simula- 
tions, determinantal |32| and world- line algorithms Q. 

A. Determinantal Monte Carlo 

For the determinantal algorithm one separates kinetic 
and potential energy terms in the Hamiltonian into the 
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product of two exponentials, and decouples the interac- 
tion term by a discrete Hubbard-Stratonovich transfor- 
mation introducing auxiliary Ising variables which we call 
fjLj(i) here @: 

1 x e A ( cr j(»))Mj(»)("jT-"ji)- AT 2 J ''' (»jT+"3l) 

2 ^ 

«w^A(tr,-(i)) = e Ar ^« (36b) 

Uj(i) ee U{<jj{i)) = U- 2guj a j (i) (36c) 

In contrast to the ordinary Hubbard model, the param- 
eter A here is not constant but depends on the local a 
variable. For the transformation Eq. (36) to be valid it 
is necessary that Uj(i) is positive for all values of <Jj(i), 
i.e. U > 2gujQ. Next one takes the trace over fermion 
degrees of freedom analytically to obtain the fermion de- 
terminant, and the Monte Carlo simulation proceeds by 
sampling the Ising spin degrees of freedom a and /i at 
each space-time site. Even though negative weights may 
occur in this formulation, we do not expect that they will 
be very significant in the dilute regime of interest for this 
model. 

The path integral formulation provided by Eq. (35), 
(36) makes the nature of this dynamic Hubbard model 
particularly apparent. The Hubbard U here has space- 
time fluctuations, with possible values Uj(i) = U — 2gu>Q 
and Uj(i) = U + 2gu!Q. This fluctuating U corresponds 
to the different values that the on-site Coulomb repul- 
sion between 2 electrons will take depending on the rel- 
ative state of these electrons, and embodies the physics 
of intra-atomic electronic correlation (at least for non- 
degenerate atomic orbitals). In a more realistic dynamic 
Hubbard model the Hubbard U will take a continuum 
of different values. The energy scale that determines the 
fluctuations in U, luq in this case, is a one-electron energy 
scale that reflects the cost in one-electron energy as the 
electrons sample the various atomic states to reduce the 
magnitude of their intra-atomic Coulomb repulsion. 

The determinantal algorithm can be used for the lattice 
problem as well as for impurity problems, and as part of 
the dynamical mean field theory solution of the model. 
This will be deferred to future work. Here we will instead 
use the world line qmc algorithm. 

B. World line Monte Carlo 

The world line Monte Carlo algorithm can be used if 
the system is one-dimensional so that no negative- weight 
problems arise. The partition function is written as 

Z = Tr f[ e~ At ^ Hj e- ArH ^ e^"^ (37) 



where the kinetic energy part of the Hamiltonian was de- 
coupled in terms involving even and odd sublattices. The 
trace in Eq. (37) is performed by introducing intermedi- 
ate states in the spin a z representation and the fermion 
occupation number representation in the usual way. In 
addition to moving fermion world lines and flipping indi- 
vidual a spins we also use composite moves consisting in 
moving a fermion world line and flipping the spins at the 
sites where the fermion occupation is changing. These 
moves are necessary to achieve equilibration in the strong 
coupling regime. 

Figure 9 shows typical world line configurations for two 
holes and the associated boson field in a strong coupling 
regime, with g = 3 and ojq — 2. We start the holes far 
from each other, in the first snapshot shown (a) after 
several hundred sweeps they are still far apart and the 
world lines are rather straight, corresponding to large 
hole effective mass. After several more hundred sweeps 
the holes bind in a bipolaron, as seen in Figure 9b. The 
bipolaron has a smaller effective mass, as indicated by the 
larger transverse motion of the world lines in the time- like 
direction. These pictures clearly show that upon pairing 
the carriers become more mobile in this model. 

The relation between pair formation and increase in 
pair mobility is shown even more clearly in Fig. 10. Fig. 
10 (a) shows the kinetic energy of a pair as a function 
of Monte Carlo sweeps. Each Monte Carlo 'step' in this 
figure gives an average over 30 consecutive sweeps. It can 
be seen that after approximately 100 steps the kinetic 
energy becomes lower. At the same time, as Fig. 10(b) 
shows, the average distance between the holes decreases 
dramatically as the pair is formed. 

When the number of holes is increased in the system it 
is found in the BCS solution in the antiadiabatic limit 
that the tendency to pairing decreases, the coherence 
length of the pairs increases until they eventually dis- 
sociate at a critical hole concentration. Similar behavior 
is found in exact diagonalization of finite systems in the 
antiadiabatic limit. We find here that similar behavior 
is seen qualitatively in Monte Carlo simulations for finite 
frequency. As an example, Fig. 11 shows snapshots of 
configurations for 6 holes in a 20-site system, i.e. hole 
concentration n./, — 0.15. The system is started in a dis- 
ordered configuration, after several thousand sweeps it 
could be seen in snapshots such as Fig. 11a that three 
well-defined pairs are formed; however, continuing the 
run, configurations like Fig. lib appear, where pairs 
overlap and the distance between members of a pair (i.e. 
the coherence length) increases. Continuing this run the 
pairs dissociate completely, later they form again. These 
snapshots suggest that for these parameters the system 
is close to the pair unbinding transition (it is not clear 
on which side) . The dependence of critical hole concen- 
tration on frequency ujq is unknown and an interesting 
subject for further study. The fact that the effective at- 
traction increases for finite frequency in the exact diag- 
onalization study suggests that the critical hole concen- 
tration may be larger for finite frequency than in the 
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cl>o — > °o limit. 

To detect a superconducting transition in Monte Carlo 
simulations it is considerably simpler to use a grand 
canonical ensemble formulation as in the determinantal 
Monte Carlo methods; in world line Monte Carlo, mea- 
surement of pairing correlation functions would involve 
breaking world lines which leads to large fluctations |33| . 
We can however get some information on pair binding 
with the world line method by consideration of the kinetic 
energy. As seen in figure 8, the kinetic energy increases 
gradually in the 4-site system as the on-site repulsion in- 
creases. Figure 12 shows the behavior of kinetic energy 
from Monte Carlo simulations on lattices of size N = 8 
and N = 12, as well as for the static (conventional) Hub- 
bard model (g = 0) . For the static Hubbard model there 
is a small dependence of kinetic energy on U, which is 
due to finite size effects. Instead, for the dynamic Hub- 
bard model there is a large increase in kinetic energy 
as U increases from small values, due to the progressive 
unbinding of the pair. For sufficiently large U the pair 
unbinds and the dependence of kinetic energy on U is 
weak as in the static Hubbbard model. For the A = 12 
cases a fairly sharp kink in the kinetic energy indicates 
the transition point. 

Similarly the existence of pairing can be seen in 
density-density correlation functions. Figure 13 (a) and 
(b) show on-site and nearest neighbor hole-hole density 
correlations for the static Hubbard model and the dy- 
namic cases of Figure 12. The on-site correlation is much 
larger in the dynamic case, and approaches the static case 
values only for large U. The nearest neighbor correlation 
in the dynamic case is much larger than in the static case 
and first increases as U increases, due to the rapid de- 
crease of the on-site correlation. Note also that for large 
U the nearest neighbor correlation is still considerably 
larger than in the static model, indicating that when the 
on-site double occupation is essentially suppressed, retar- 
dation gives rise to an effective nearest neighbor attrac- 
tion. This can be easily understood from second order 
strong coupling perturbation theory. 

Electrons behave very different from holes in this 
model. In figure 14 we show snpshots of hole worldlinc 
configurations when the band is almost full with holes, 
i.e. almost empty with electrons. This is the mirror 
image of the case shown in Fig. 9. Here we start the 
simulation with two electrons on the same site, and af- 
ter some sweeps the electrons separate. Furthermore, in 
contrast to Fig. 9, the quasiparticle world lines show 
much larger fluctations in the time direction indicating 
the smaller electron effective mass. The contrast between 
Fig. 14 and Fig. 9 clearly displays the intrinsic electron- 
hole asymmetry of this dynamic Hubbard model. Kinetic 
energy in this case shows almost no dependence on [/, as 
expected, in contrast to the case shown in Fig. 12. 

In summary, the results of these Monte Carlo simula- 
tions support the picture obtained from exact diagonal- 
ization of small systems: the dynamic Hubbard model is 
an effective way to obtain hole pairing driven by kinetic 



energy lowering in repulsive fermion systems. 

V. COMPARISON WITH RESULTS FOR A 
HOLSTEIN-LIKE MODEL 

The conventional electron-boson models studied in the 
past involve coupling of a boson degree of freedom to the 
electronic charge density rather than to the double occu- 
pancy. Even though it doesn't necessarily follow, in their 
simplest form these models are electron-hole symmetric. 
We consider here one such model with site Hamiltonian 

H t = LO a l x + guj Q al[n^ + n il - 1] + Un^n,^ (38) 

as a generic model in that class. This model should be 
similar to the Holstein model p5|| , where the spin- 1/2 
degree of freedom is replaced by a harmonic oscillator. 
Diagonalization of the site Hamiltonian yields eigenvalues 

e(0) = e(2) = -e(0) = -e(2) = -cj oy /l + g 2 (39a) 

e(l) = -e(l) = -u (39b) 
and eigenvectors of the form Eq. (4), with 

u = -UA+ ; 9 ( 40a ) 

V2\J a/TTF 

v = -^ p .ji~- - 9 (40b) 

x/2y ^tt7 

and effective on-site interaction 

[/ e// = f/-2w (Vl+5 2 -l). (41) 

The overlap matrix element between ground state wave 
functions is 

(42) 

and ranges between 1 for g = to l/v2 for g — > oo, so 
that it never becomes small as in the previous case. 

Figure 15 shows results of exact diagonalization for an 
A = 4 cluster. As a function of freqency, the effective 
interaction becomes less attractive as ujq decreases, in 
contrast to the behavior found for the dynamic Hubbard 
model. Furthermore, the kinetic energy increases as loq 
increases and the effective attraction increases, as shown 
in Fig. 15b. Hence in this model pairing gives rise to ki- 
netic energy increase, and the pair condensation energy 
originates in the larger potential energy decrease, which 
is precisely opposite to the behavior in the dynamic Hub- 
bard model. 
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Figure 16 shows typical world line configurations for 
two holes in this model. Similarly to Fig. 9, we start the 
holes far from each other, in the first snapshot they are 
still separate and after several sweeps a pair is formed. 
Here the world lines for the pair are rather straight, cor- 
responding to large effective mass, while the single holes 
exhibit larger transverse motion of the world lines indi- 
cating lighter quasiparticles. This is qualitatively differ- 
ent to the behavior in the dynamic Hubbard model (fig. 
9), where the carriers became lighter when they paired. 

The relation between pair formation and pair mobility 
is also shown in Figure 17, to be compared with Figure 10 
for the dynamic Hubbard model. Initially the hole-hole 
distance is large and the kinetic energy is low. When the 
pairs form after approximately 100 Monte Carlo 'steps' 
the hole-hole distance decreases drastically and the ki- 
netic energy increases, again in qualitative contrast with 
the behavior found in the dynamic Hubbard model. 

Finally, Figure 18 shows the behavior of the average 
kinetic energy and of density-density correlations versus 
on-site repulsion U for the Holstein-like model. As in- 
dicated by Fig. 15, the effective interaction becomes at- 
tractive for g = 3 and cjq = 1 when U *~ 4. This is 
confirmed by the results of Fig. 18. The on-site density- 
density correlation increases sharply as U is decreased 
below 4, indicating pair formation. At the same time, 
the kinetic energy increases sharply when the pair forms, 
again in stark contrast to the behavior of the dynamic 
Hubbard model seen in Fig. 12 and 13. The nearest- 
neighbor density-density correlation first increases as U 
increases from 0, indicating that the pair wavefunction 
evolves from describing on-site pairing to more extended 
pairing, and then decreases for larger U as the pair dis- 
sociates. 

In summary, these results suggest that the conven- 
tional electron-hole symmetric models and dynamic Hub- 
bard models define two rather different 'universality 
classes'. Both types of models can describe pairing, with 
qualitatively different features. 

VI. RELATION TO REAL ATOM 

In the dynamic Hubbard model considered in this pa- 
per, the on-site repulsion takes the values U — 2gojQ and 
U + 2gujQ when the auxiliary spin at the site points up 
and down respectively in a a z representation. More gen- 
erally, for the spin in a superposition of these states the 
on-site repulsion will take values intermediate between 
these extremes. The 'effective' on-site repulsion defined 
by 

U e ff(site) = £7(2) + £7(0) - 2£7(1) (43) 

with £7(n) the site energy with n electrons (or holes) is 
simply U. The reason a fluctuating U is needed to repre- 
sent a real atom is that the wavefunction of 2 electrons in 
an orbital is not simply the product of the single-electron 



wavefunctions in the singly occupied atom |3J , but rather 
a superposition 

#(ri,r 2 ) = J2 Cnm ¥n{ri)(p m {r 2 ) (44) 

where {</?n(<")} is a complete set of single-electron wave 
functions. The fluctuating values of U can be thought 
of as the different values that the electron-electron re- 
pulsion will take for one electron in <p n (f) and the other 
electron in (p m (r), for all n, m for which C nm is not 
zero. The frequency loq represents the energy scale of 
electronic excitations in the atom, i.e. the eigenenergies 
of the wavefunctions {ip n (r)} . 

More specifically, for the particular case of Is orbitals 
in a hydrogenic atom of ionic charge Z, the 'bare' on-site 
repulsion for two electrons in the Is orbital is 

Ui = 17ZeV (45) 

This corresponds in our model to the on-site repulsion 
when the boson is not allowed to relax, Eq. (8), or ap- 
proximately U + 2gujQ. In the Hartree approximation, 
the orbital expands to Z = Z — 5/16 upon double occu- 
pation, and the repulsion between two electrons in these 
expanded orbitals is Eq. (46) with Z replaced by Z, i.e. 

U 2 = 17 (Z - —)eV = U l - b.ZleV (46) 
16 

This would roughly correspond to the 'minimum' on-site 
repulsion in our model, U — 2gujQ. Finally, the effective 
on-site repulsion in the Hartree approximation, taking 
into account the cost in single-particle energy upon or- 
bital expansion, is 

U 3 = 17(Z - —)eV = Ui- 2.66eV (47) 

precisely halfway between the values Eq. (46) and (47), 
and this would correspond to the effective site U in our 
model, which is just U , also precisely halfway between 
the minimum and maximum values. Summarizing, 

Ui = 17ZeV =U + 2gLU (48a) 

U 2 = Ui- 5.3 W = U- 2guj {) (48b) 

U 3 = Ui- 2.66eV = U (48c) 
so that for this case we have simply 

gu = 1.33eV (49) 

independent of Z. 

The frequency ojq is related to excitation energies of the 
atom, hence we expect the dependence on ionic charge 

uj = cZ 2 (50) 







which implies from Eq. (50) that g increases as the ionic 
charge Z decreases. This is in accordance with the fact 
that the overlap 



S =< Oil >= 



1 



(51) 



decreases as g increases; in the Hartree approximation, 
this overlap corresponds to the overlap of the expanded 
single electron orbital with the non-expanded one and is 
given by 



S = 



(1 



(52) 



32Z ' 



which decreases to zero as Z — > 0.3125. 

Strictly speaking our dynamic Hubbard model will be 
a valid representation of the real atom only in the pa- 
rameter regime where 



U - 2goj Q > 



(53) 



because the atomic Coulomb integral for any two orbitals 
(p n , ip m has to be positive. Even with the constraint Eq. 
(54) a wide range of parametes in the model exists where 
pairing will occur, as can be inferred from the numerical 
results in the previous sections. For example, from Fig. 
la we see that for lj = 0.5, U = 8, pairing occurs for 
g > 4; the condition Eq. (54) is satisfied in this case up to 
.9 = 8, and for g ~ 4 the fluctuations in the on-site U are 
about 50%. As the frequency gets smaller, the relative 
fluctuations in U needed to obtain pairing decrease. For 
example, from Figure 2a we find that for U = 4, g = 
4, pairing occurs for loq > 0.12, which corresponds to 
fluctuations in U of only 25% (between [7 = 3 and U — 
5). If the U was not fluctuating but fixed, no pairing 
occurs in the model unless U < 0. In other words, the 
'equivalent V in a model with fixed U is not only smaller 
than the average U in the fluctuating case but it is smaller 
than the smallest value that the fluctuating U attains in 
these cases. 

Note that as Z in the atom decreases, the on-site bare 
U decreases (Eq. (46)), the parameter ojq should decrease 
according to Eq. (51), and correspondingly g should in- 
crease (Eq. (50). As seen in Figure 2b, as U decreases 
a smaller loq is required to give pairing, and as seen in 
Figure 2a as g increases also a smaller ujq is needed for 
pairing. We conclude from our results for the model sys- 
tem and the relationship with the real atom that smaller 
values of Z yield the most favorable conditions for pairing 
in this dynamic Hubbard model. 



VII. DISCUSSION 

We have studied numerically some properties of a dy- 
namic Hubbard model, where the value of the on-site 
repulsion U depends on the state of an auxiliary boson 



degree of freedom. In the model studied in this paper 
the boson is a spin- 1/2 degree of freedom, with excita- 
tion energy loq. It will be of interest to study other similar 
models with other boson degrees of freedom such as other 
versions of the spin 1/2 model JjJ], higher spin variables 
or harmonic oscillators, or purely electronic models with 
more than one orbital per site we believe the qual- 
itative physics will be similar. Furthermore, 'extended' 
dynamic Hubbard models, with more than on-site inter- 
actions, should be interesting to study. 

These dynamic Hubbard models map onto the Hub- 
bard model with correlated hopping in the antiadiabatic 
limit luq — > oo, which is known to lead to pairing of holes 
and superconductivity for sufficiently large coupling con- 
stant g ; the purpose of this paper was to determine 
whether pairing still exists for finite ujq. Furthermore in 
the antiadiabatic limit pairing is known to occur through 
kinetic energy lowering, and we examined whether the 
same physics occurs for finite u>q. 

From both the exact diagonalization and the quantum 
Monte Carlo results we concluded that the same physics 
of the antiadiabatic limit persists for finite, and even 
small, loq. Furthermore, the parameter regime where 
pairing occurs is larger for small ujq than for loq — > oo. 
When ujq is small the kinetic energy is much lower than 
in the antiadiabatic limit, yet the magnitude of kinetic 
energy lowering upon pairing is similar to that in the 
antiadiabatic limit. 

The ujq — ► oo limit of the model (Hubbard model with 
correlated hopping) is useful because its physics is rather 
transparent and because it allows for much simpler ana- 
lytic and numerical treatments. However, strictly speak- 
ing the dynamic Hubbard model considered here is only 
a realistic representation of a real system for parameters 
where U — 2gtuo > 0, which certainly does not hold in 
the antiadiabatic limit. Hence it is essential to estab- 
lish that the properties of the model for small ujq and for 
luq oo are similar if one is to use the results obtained 
from the u>o — > oo limit to understand the properties of 
a real system. We found that the model can give rise 
to pairing even in parameter ranges where the fluctuat- 
ing U attains only positive values, which implies that an 
'equivalent' fixed U in those cases would be smaller than 
the lower bound of the range within which U fluctuates. 

We contrasted the behavior found in the dynamic 
Hubbard model with that of an electron-hole symmet- 
ric Holstein-like model. In the latter model, which we 
suggest is representative of a wide range of model Hamil- 
tonians that have been considered in the past to describe 
superconductivity, the physics found is qualitatively dif- 
ferent: pairing is associated with lowering of potential 
energy and increase in kinetic energy, opposite to the 
behavior found in the dynamic Hubbard model. We 
suggest that these two models, each representative of 
an entire class of model Hamiltonians, are two different 
paradigms by which superconductivity can be achieved 
p6| . Whether either or both occur in nature is an un- 
settled question. The theory of hole superconductivity 
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proposes that only the paradigm represented by the dy- 
namic Hubbard model occurs in real materials [^flp] . 

We also discussed briefly the relation between the 
dynamic Hubbard model considered here and a real 
hydrogen-like atom for Is electrons in the Hartree ap- 
proximation. Clearly such relation should be qualita- 
tively similar for more accurate representations of the 
two-electron wavefunction such as the Hylleraas wave 
function as well as for electrons in other atomic or- 
bitals. We found that a smaller value of the ionic charge 
Z yields more favorable conditions for pairing for sev- 
eral different reasons: 1) it leads to smaller on-site re- 
pulsion U, 2) it leads to larger coupling g, which leads 
to larger 'dressing' of quasiparticles in the normal state, 
and to larger 'undressing', hence larger energy lowering, 
as quasiparticles pair, and 3) it leads to smaller frequency 
scale loq, which according to the results of this paper is 
favorable to pairing as U becomes smaller and g becomes 
larger. In addition, smaller Z is also favorable because 
it leads to larger orbital overlap between atoms, hence 
larger bare hopping t, which increases the overall scale of 
the pairing interaction calculated in this paper plj ). It 
will be interesting to perform detailed analysis of the con- 
nection between this and other dynamic Hubbard models 
and electrons in various orbitals in real atoms. 

The results presented here are only a first step in the 
understanding of this and other dynamic Hubbard mod- 
els. Analytically, both strong and weak coupling expan- 
sions should be feasible and of interest. Powerful numer- 
ical techniques that have been extensively used for the 
static Hubbard model and other related models such as 
the Holstein model can and should be brought to bear on 
this class of models. In particular, density matrix renor- 
malization group |35| ] and determinantal Monte Carlo 
methods should allow for the study of larger systems as 
a function of model parameters and hole concentration 
to determine the range of parameters where hole pair- 
ing occurs. The dynamical mean field method combined 
with a Monte Carlo 'impurity' method should be a very 
fruitful approach to deal with this class of models |S7| . 
In particular, it will be of great interest to understand 
quantitatively the processes of spectral weight transfer 
in one- and two-particle Green's functions that are ex- 
pected to occur in this class of models upon transition to 
the superconducting state |3^,^l|, which are of interest 
in connection with photoemission p8| and optical exper- 
iments [p9| in superconducting materials. 
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FIG. 1. Effective interaction U e ff for N = 4 cluster versus 
coupling constant g and various values of the on-site repulsion 
U for (a) u>o = 0.5 and (b) loo — 2 (full lines). The dashed 
lines and dotted lines give the results in the luo — » oo limit for 
the N — 4 cluster and for the infinite chain respectively. For 
fixed g, increasing U corresponds to increasing value of U e ff- 



FIG. 2. Dependence of U e f / on uiq for N 
dashed lines give the limiting values u>o — *• oc 



4 cluster. The 



FIG. 3. Phase diagram for (a) N = 2 and (b) N = 4 clus- 
ters. In the region labeled NON-SC, U e ff > for all values 
of u>o', in the region labeled SC, a range of uio exists where 
Ueff < 0. Below the solid line, U e ff < in the antiadiabatic 
limit luq — > oo. 



FIG. 7. Difference between kinetic energy of a pair and 
kinetic energy of two holes in the 4-site chain (a) versus g for 
fixed wo and (b) versus luo for fixed g. The dashed and dotted 
lines in (a) give the results in the luo — > oo limit for the N = 4 
cluster and the infinite chain respectively; the dashed lines 
in (b) give the results in the luo — » oo limit for the N = 4 
cluster, (c) Ratio of kinetic energy lowering to pair binding 
energy in 4-site chain versus g for (7 = 2 and various values of 
u>o- The dashed line gives the results for u>o = oo. Note that 
the kinetic energy lowering upon pairing for small frequency 
is considerably larger for finite luo than in the loo = oo limit. 



FIG. 8. Single particle (full lines) and pair (dashed lines) 
kinetic energy versus U for g = 3 and g = 4 for (a) luo = 0.5 
and (b) loo = 2. The results for luo — > oo are also shown. 
Both the single particle and the pair kinetic energies are sub- 
stantially lower for small luo than for luq — oo. The transition 
points where U e ff changes sign for g — 3 and g = 4 are re- 
spectively U — 3.6 and U = 3.77 for u>o = oo, U = 6.91 and 
U = 10.1 for lu = 2 and U = 5.57 and U = 7.77 for lu = 0.5. 



FIG. 9. Snapshots of Monte Carlo configurations for an 
N = 20-site lattice with L = 40 time slices. At = 0.25, 
g = 3, u>o — 2. The left-side panel show the hole world lines, 
the right-hand panel indicates the boson (spin) configuration. 
At the sites where the hole occupation is 1 (2), the boson 
configuration is denoted by p (6) if the boson state is |— >, 
which is the low energy configuration, and by n if it is + >. 
In (a), after a few hundred sweeps, the holes are separate 
and heavy, in (b), after several more hundred sweeps, the 
holes are bound and lighter (world lines fluctuate more in 
time direction). 



FIG. 4. Optimal frequency luo that gives rise to pairing at 
the phase boundaries (dashed lines) of figure 3. As g decreases 
too increases and the phase boundary approaches the one in 
the antiadiabatic limit. 

FIG. 5. Results for kinetic energy in the infinite chain in 
the antiadiabatic limit, (a) Kinetic energy of two unbound 
holes (full line) and of a hole pair (dashed lines) versus g for 
various values of U. As g decreases, the dashed line joins the 
full line (as indicated by the symbols) when pairs unbind, at 
g = g c . g c is 0.5, 1 and 2 for U = 0.8, 2 and 3.2 respectively, 
(b) Pair binding energy (full lines) and kinetic energy lowering 
(dashed lines) for the infinite chain versus g for various U. At 
g c , both quantities go to zero, (c) Ratio of kinetic energy 
lowering to pair binding energy versus g. 

FIG. 6. Comparison of results for kinetic energy of iV = 4 
chain and infinite chain in the antiadiabatic limit. In (a), the 
kinetic energy of 2 holes in the N = 4 cluster (dash-dotted 
line) joins the infinite chain kinetic energy of the bound pair 
for large g (dashed line) and the kinetic energy of unbound 
holes for small g (full line), (b) shows that finite size effects 
similarly enhance the magnitude of pair binding energy and 
of kinetic energy lowering. 



FIG. 10. Kinetic energy of two holes (a) and hole-hole dis- 
tance (b) as function of 'Monte Carlo time'. The unit of time 
is 30 Monte Carlo sweeps. N = 20, L = 40, At = 0.25, g = 3, 
uio = 1. After about 100 steps the hole- hole distance decreases 
drastically as the holes become bound (b), and at the same 
time the kinetic energy becomes lower (a). 

FIG. 11. Snapshots for 6 holes in a 20-site chain. L = 40, 
At = 0.25, g = 3, loo = 2. In the first picture (a) three 
well-defined pairs are seen, in (b) the pairs overlap and the 
distance between members of a pair ('coherence length') in- 
creases. 

FIG. 12. Kinetic energy versus on-site repulsion U from 
Monte Carlo simulations of dynamic Hubbard model (a), (b), 
and of static (conventional) Hubbard model (g — 0) (c). 

FIG. 13. (a) On-site and (b) nearest neighbor den- 
sity-density correlations for dynamic and static Hubbard 
model, TV = 12 sites. 
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FIG. 14. Snapshots of Monte Carlo configurations when 
the band is almost full with holes, with only two electrons. 
The convention and parameters are the same as in Fig. 9. The 
simulation is started with the two electrons on the same site 
(top panel), after a few sweeps (bottom panel) the electrons 
separate and are light, as illustrated by the large fluctuations 
of the world lines in the time direction. 

FIG. 15. Exact diagonalization results for Holstein-like 
model Eq. (39). (a) Effective interaction and (b) difference 
between kinetic energy of a pair and kinetic energy of two 
holes in the 4-site chain versus frequency wo for U — 4 and 
various values of g. The kinetic energy increases rapidly as 
the effective interaction becomes more attractive, in contrast 
to the behavior found in the dynamic Hubbard model. 

FIG. 16. Snapshots of Monte Carlo configurations for the 
Holstein-like model for an N — 20-site lattice with L — 40 
time slices. At = 0.25, g = 3, a>o = 2, U = 4. Same conven- 
tions as in Fig. 9. After a few initial sweeps the holes are still 
separate and light (upper panels), after several more sweeps 
(lower panels) the holes are bound and heavier (world lines 
fluctuate less in time direction), in contrast to the behavior 
seen in Fig. 9. 

FIG. 17. Kinetic energy of two holes (a) and hole-hole 
distance (b) as function of 'Monte Carlo time' for the Hol- 
stein-like model. The unit of time is 15 Monte Carlo sweeps. 
N = 20, L = 40, At = 0.25, g = 3, w = 1, U = 3. After 
about 100 steps the hole-hole distance decreases drastically as 
the holes become bound (b), and at the same time the kinetic 
energy increases (a), in contrast to the behavior seen in Fig. 
10. 

FIG. 18. (a) Kinetic energy and (b) on-site and nearest 
neighbor density- density correlation versus on-site repulsion 
U from Monte Carlo simulations of the Holstein-like model 
with N = 12. In (a) results for N = 8 are also shown, for (b) 
the size dependence is negligible. As U decreases and pairs 
form, the on-site correlations increase and the kinetic energy 
increases, in contrast to the behavior seen in Figs. 12 and 13 
for the dynamic Hubbard model. 



13 



I I I I I I I I I I I I I I I I I I I I I I I I I I I 



(a) U=4 





Figure 



12 
10 
8 
6 
4 
2 




j i i i 



i i i i 



i i i i 



i i i i 



(a) N=2 



MM 

/ =1 
/- 



/ 



NON-SC / - 

/ 

/ _ 

sc □ 



/ 



/ 




I I I I I I I 



I I I I 



I I I I 



I I I I' 



1 2 3 4 5 



15 



10 



MM 



MM 



MM 



MM 



(b) N=4 



MM 



/ - 
/ - 



/ 



NON-SC 



/ 



SC 



l-H I I I I I I I I I I I I I I I I I I I I I I 

1 2 3 4 5 



Figure 3 



0.0 
■0.5 
1.0 
■1.5 
•2.0 
•2.5 
•3.0 



1 II J 
- n=4 


i i i i 


II II Ml II. 




^ = ~V A 

- ^ ^bntiad. - 
~ ^ z 

z (a) co =0.5 : 


- q-4 




- ^ 
-\ i i i 


c 

i i i i 


- " -- 

single part. - 

- -pairl 
i r i i i i i i i 







4 
U 



8 



0.0 



0.5 



1.0 



1.5 



1 1 


i i 

=4 


i i i i 


1 1 1 1 

(b) co 


i i i i 
V2 ~_ 


~ q= 


=3 




— antiad. - 






^ , 




- / 


/ 






g=4 -_ 


y 




y 


g=3 ^ 




/ — single part. - 








- -pair 




i i 


i i 


1 1 1 1 


i i i i 


i i i i 







4 
U 



8 



A 
I— 

V 



I I I I I I I I I I I I I I I I I I I 



I I I I 



(a) Kinetic 
energy 




i i i i 



i i i i 



i i i i 



50 100 150 200 
MC time 



A 
V 



15.0 



i i i i i i i i i i i i i i i i i i i 



(b) hole-hole 
distance 




50 100 150 200 
MC time 



Figure 10 



********************** ********************** 

* * 

*+++++++b++p+p+++++b+* 

* * 

*+++++++b++p+p+++++b+* 

* * 

*+++++++pp+p+p+++++b+* 

* * 

*-++++++pp+p+p+++++b+* 

* * 

*+++++++pp+b+++++++b+* 

* * 

*+++++++pp++b++++pp++* 

* * 

*+++++++pp+n+p+++pp++* 

* * 

*+++++++pp++pp+++pp++* 

* * 

*+++++++pp+-+b++++b++* 

* * 

*+++++++pp+++pp++pp++* 

* * 

*+++++++pp+++b+++pp++* 

* * 

*+++++++b++++b+++pp++* 

* * 

*++++++pp++++b+++pp++* 



* 


# 


- + 


# 


* 


* 


# 


- + 


# 


* 


* 


# 


- + 


# 


* 


* 


# 


- + 


# 


* 


* 


+- 


- + 


# 


* 


* 


+- 


- + 


# 


* 


* 


-+ 


- + 


# 


* 


* 


-+ 


- + 


# 


* 


* 


-+ 


-+ 


# 


* 


* 


-+ 


# 


# 


* 


* 


-+ 


# 


# 


* 


* 


-+ 


# 


-+ 


* 


* 


-+ 


+- 


-+ 


* 


* 


-+ 


+ - 


-+ 


* 


* 


-+ 


+ - 


-+ 


* 


* 


-+ 


+- 


-+ 


* 


* 


-+ 


# 


-+ 


* 


* 


-+ 


# 


# 


* 


* 


-+ 


# 


# 


* 


* 


-+ 


-+ 


-+ 


* 


* 


-+ 


-+ 


-+ 


* 


* 


-+ 


# 


-+ 


* 


* 


# 


# 


-+ 


* 


* 


# 


# 


-+ 


* 


* 


# 


# 


-+ 


* 


* 


-+ 


# 


-+ 


* 



********************** ********************** 

* * 

*++++b+pp++++++++++b+* 

* * 

*++++b+pp++++++++++b+* 

* * 

*++++pppp+++++++++np+* 

* * 

*++++b+pp+++++++++p+p* 

* * 

*++++pppp+++++++++p+p* 

* * 

*++++pppp+++++++++pn-* 

* * 

*++++b+pp++++++++pp++* 

* * 

*++++b+pp+++++++++b++* 

* * 

*++p+pp+p+++++++++b++* 

* * 

*++p+pp+p+++++++++np+* 

* * 

*++p+b++p+++++++++np+* 

* * 

*++ppp++p++++++++++b+* 

* * 

*++ppp++p++++++++++b+* 



* 


# -+ 


# * 


* 


# -+ 


# * 


* 


# -+ 


# * 


* 


# +- 


# * 


* 


# +- 


+- * 


* 


+-+- 


+- * 


* 


+-+- 


-+ * 


* 


# -+ 


- +* 


* 


# -+ 


- +* 


* 


+ — + 


- +* 


* 


+ — + 


- +* 


* 


+ — + 


-+ * 


* 


+ — + 


# * 


* 


# -+ 


+- * 


* 


# -+ 


+- * 


* 


# -+ 


# * 


* 


-+- + 


# * 


* 


- +- + 


# * 


* 


- +- + 


+- * 


* 


- -+ + 


+- * 


* 


- -+ + 


+- * 


* 


- # + 


+- * 


* 


-+- + 


# * 


* 


+ — + 


# * 


* 


+ — + 


# * 


* 


+ — + 


# * 



Figure 11 



•1.8 
■2.0 
•2.2 
•2.4 



_i i i i I i i i i I i i i i I i i i i. 
I (a) qh, co =lo.5 1 




I — I I I I I I I I I I I I I I I I I I I — I 




•3.2 



.1 i i i I i i i i I i i i i I i i i i. 



•3.4 |t- (c) q-0 
■3.6 
3.8 

•4.0 I — 



N-12 

I— N=8 —I 



£P x X X X — j( 



.4 2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 

" 2 4 6 8 
U 



Figure 12 




Mill 



0.06 



m 



0.04 



0.02 



0.00 



i i i i 



i i i i 



I I I H 



C] t (b) <n 1 n 1+1 >_ 



$ I A ^ 



J>q=3 9 co =0.5 
-©q=4 9 co =0.5 * $ ^~ 



K Xq=0 



i i i i 







4 
U 



8 



Figure 13 



********************** 

*######### ##########* 
*######### ##########* 
*######### ##########* 
*######### ##########* 
*######### ##########* 
*######### ##########* 



********************** 

* * 

*bbbbbbbbb+bbbbbbbbbb* 

* * 

*bbbbbbbbb+bbbbbbbbbb* 

* * 

*bbbbbbbbb+bbbbbbbbbb* 



*######### 
*######### 
*######### 
*######### 
*######### 
*######### 
*######### 



*######### 
*######### 
*######### 
*######### 
*######### 
*######### 
*######### 
*######### 
*######### 
*######### 
*######### 
*######### 



##########* 
##########* 
##########* 
##########* 
##########* 
##########* 
##########* 
##########* 
##########* 
##########* 
##########* 
##########* 
##########* 
##########* 
##########* 
##########* 
##########* 
##########* 
##########* 
##########* 
##########* 
##########* 



*bbbbbbbbb+bbbbbbbbbb* 

* * 

*bbbbbbbbb+bbbbbbbbbb* 

* * 

*bbbbbbbbb+bbbbbbbbbb* 

* * 

*bbbbbbbbb+bbbbbbbbbb* 

* * 

*bbbbbbbbb+bbbbbbbbbb* 

* * 

*bbbbbbbbb+bbbbbbbbbb* 

* * 

*bbbbbbbbb+bbbbbbbbbb* 

* * 

*bbbbbbbbb+bbbbbbbbbb* 

* * 

*bbbbbbbbb+bbbbbbbbbb* 

* * 

*bbbbbbbbb+bbbbbbbbbb* 

* * 

*bbbbbbbbb+bbbbbbbbbb* 

* * 



********************** 
*#+#############-####* 

*##+############-####* 

*##+############-####* 
*###+##########-#####* 
*####+########-######* 
*####+########-######* 
*####+#########-#####* 
*####+#########-#####* 
*####+#########-#####* 
*#####+########-#####* 
*#####+########-#####* 
*#####+########-#####* 
*#####+########-#####* 
*####+##########-####* 
*###+###########-####* 
*##+############-####* 
*##+############-####* 
*##+############-####* 
*##+############-####* 

*##+#############-###* 
*##+##############-##* 
*##+##############-##* 
*###+#############-##* 
*###+#############-##* 
*###+#############-##* 
*###+#############-##* 



********************** 

* * 

*bpbbbbbbbbbbbbbpbbbb* 

* * 

*bbpbbbbbbbbbbbbpbbbb* 

* * 

*bbbpbbbbbbbbbbpbbbbb* 

* * 

*bbbbpbbbbbbbnpbbbbbb* 

* * 

*bbbbpbbbbbbbnbpbbbbb* 

* * 

*bbbbbpbbbbbbnbpbbnbb* 

* * 

*bbbbbpbbbbbbbbpbbnbb* 

* * 

*bbbbpbbbnbbbbbbpbbbb* 

* * 

*bbpbbbbbnbbbbbbpbbbb* 

* * 

*bbpbbbbbnbnbbbbpbbbb* 

* * 

*bbpbbbbbnbbbbbbpbbbb* 

* * 

*bbpbbbbbnbbbbbbbbpbb* 

* * 

*bbbpbbbbbbbbbnbbbpbb* 

* * 

*bbbpbbbbbbbbbnbbbpbb* 



Figure 14 



********************** ********************** 

* _ + * * * 

* _ + * * * 

* _ + * * + +++ + +++ ++ + n + + + + ++ n +* 

* _ + * * * 

* + * *++++++++++++n+++++n+* 

* _ + * * * 

* + * *++++++++++++n+++n+++* 

* _ + * * * 

* - + * *++++++++++++n+++p+++* 

* _ + * * * 

* - + * *+++++++++++n+++++n++* 

* _+** * 

* - + * *++-++++++++n++++n+++* 

* _-(-** * 

* - + * *+++++++++++n++-n++++* 

* _ + * * * 

* + * *+++++++++++p+++++n++* 

* _+** * 

* _ + * * * 

* - + * *++++-+-++++++n++n+++* 

* _ + * * * 

* _ + * * * 

* - + * *+++++++++++++p+++n++* 

* _+** * 

* - + * *+++++++++++++n+++n++* 

* _+** * 

********************** ********************** 

* +_ * * * 

* +- * *++++++++++++++nn++++* 

* +_ * * * 

* +- * *++++++++++++++pn++++* 

* +_ * * * 

* # * *++++++++++++++b+++++* 

* # * * * 

* # * *++++++++++++++fo+_+++* 

* # * * * 

* # * * * 

* # * *++++++++++++++b+++++* 

* # * * * 

* # * *++++++++++++++b+++++* 

* # * * * 

* # * *++++++++++++++b+++++* 

* # * * * 

* # * *++++++++++++++b+++++* 

* # * * * 

* # * *+++++++++++++-b+++++* 

* # * * * 

* # * *++++++++++++++b++-++* 

* # * * * 

* _+ * * ++++++++++++++ p n++++ * 

* _+ * * * 

* _+ * * ++++++++++++++ p n++++ * 

* _+ * * * 

* _+ * * * 



Figure 16 



■2.0 
■2.5 
■3.0 
■3.5 
■4.0 



i i i i 



i i i i 

- (a)Kinetic 
- energy 



i i i i 



i i i i 




.4 5 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 

50 100 150 200 
MC time 



15.0 



1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 



(b) hole-hole 
distance 




ri 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 
50 100 150 200 
MC time 



A 
I— 

V 



I I I I 



N-12 




i i i i 



-X- a 



Q I I I I I I I I I I I I I I I 

2 4 6 8 
U 



LI I I I I I I I I I I I I I I I I I IJ 



0.3 



0.2 

c 



c 

V 



0.1 



0.0 



(b) 



\ 



<n 1 n 1 > 



*-* i nin J +1> l 



2 4 6 8 

u 



Figure 18 



